#include "stdafx.h"
#include "PNGFinderEstimationAlgorithm.h"

namespace SpikingNeuronLib
{
	PNGFinderEstimationAlgorithm::PNGFinderEstimationAlgorithm(List<Connection^>^ connections, int numSynapses, int numInhibitoryNeurons, int totalNeurons, int groupSize, int maxPNGLength)
		: PNGFinder(connections, groupSize)
	{
		_convergedEvents = gcnew array<VirtualEvent>(3 * numSynapses);
		_inhibitionList = gcnew array<FiringEvent>(3 * numInhibitoryNeurons * numSynapses);

		_virtualEvents = gcnew array<VirtualEvent, 2>(maxPNGLength, totalNeurons);
	}

	int PNGFinderEstimationAlgorithm::SearchCurrentTriplet(String^ tag, int anchorNeuron, array<int>^ triplet, int maxDelay, int numPolychronousGroups,
		int groupSize, double weightLimit, int maxPNGLength, int numExcitatoryNeurons, int totalNeurons,
		Dictionary<Tuple<int, int>^, Connection^>^ connectionPairDictionary,
		Dictionary<int, Dictionary<int, List<Connection^>^>^>^ groupedConnectionDictionary, int verboseLevel)
	{
		return SearchCurrentTripletEstimated(tag, anchorNeuron, triplet, maxDelay, numPolychronousGroups,
			groupSize, weightLimit, maxPNGLength, numExcitatoryNeurons, totalNeurons,
			connectionPairDictionary, groupedConnectionDictionary,
			NumPostspikes, GroupEvents, _virtualEvents, Links, _convergedEvents, _inhibitionList, verboseLevel);
	}

	int PNGFinderEstimationAlgorithm::SearchCurrentTripletEstimated(String^ tag, int anchorNeuron, array<int>^ triplet, int maxDelay, int numPolychronousGroups,
		int groupSize, double weightLimit, int maxPNGLength, int numExcitatoryNeurons, int totalNeurons,
		Dictionary<Tuple<int, int>^, Connection^>^ connectionPairDictionary,
		Dictionary<int, Dictionary<int, List<Connection^>^>^>^ groupedConnectionDictionary,
		array<int>^ numPostspikes, array<FiringEvent>^ groupEvents, array<VirtualEvent, 2>^ virtualEvents,
		List<LinkedEventPair^>^ links, array<VirtualEvent>^ convergedEvents, array<FiringEvent>^ inhibitionList, int verboseLevel)
	{
		const int minConvergentSpikes = 2;      // number of spikes needed to fire a cell
		const int latency = 4;                  // latency from the moment spikes arrive to the moment the neuron fires
		const int temporalConvergenceWidth = 2;
		const int inhibitoryConvergenceWidth = 5;
		const int slacke = 0;
		const int slacki = 0;
		const int minPNGPathLinks = 7;          // minimal length of a group

		// MG originally used unassigned variable 'timing' that was assumed to be set to the max timing for each iteration
		// based on the order of iterating through the delays
		// this new implementation cannot assume that delays are iterated in a particular order
		// therefore the new implementation explicity sets the max timing for each outer iteration of time, t
		int maxTiming = 0;

		// reinitialise
		int numFired = groupSize;
		int numInhibitionList = 0;
		int maxLinkLevel = 1;
		links->Clear();

		// fetch the triplet connections
		array<Connection^>^ tripletAnchorConnections = gcnew array<Connection^>(triplet->Length);
		for (int i = 0; i < triplet->Length; i++)
		{
			Connection^ connection = connectionPairDictionary[gcnew Tuple<int, int>(triplet[i], anchorNeuron)];
			tripletAnchorConnections[i] = connection;
		}

		// find the longest connection for the triplet
		int Dmax = 0;
		for (int i = 0; i < groupSize; i++)
		{
			Connection^ connection = tripletAnchorConnections[i];
			if (Dmax < connection->Delay)
			{
				Dmax = connection->Delay;
			}
		}

		// initialize the data structures
		// for each triplet neuron...
		for (int groupIndex = 0; groupIndex < groupSize; groupIndex++)
		{
			int neuronIndex = triplet[groupIndex];
			Connection^ anchorConnection = tripletAnchorConnections[groupIndex];

			// add the neuron (and firing time) to the new PNG
			groupEvents[groupIndex] = FiringEvent(Dmax - anchorConnection->Delay, neuronIndex, EventLabel::Foreground);

			// for each outgoing connection from neuronIndex
			Dictionary<int, List<Connection^>^>^ outgoingConnections = groupedConnectionDictionary[neuronIndex];
			for each(int delay in outgoingConnections->Keys)
			{
				for each(Connection^ connection in outgoingConnections[delay])
				{
					// if the connection is to an inhibitory neuron
					// or the connection is to a strong excitatory neuron that will fire after the anchor...
					if ((connection->PostNeuron >= numExcitatoryNeurons) || ((connection->Weight > weightLimit) && (delay >= anchorConnection->Delay)))
					{
						// tentatively add to the group
						//int timing = groupEvents[i].Time + delay + 1;     // in the original, delay is zero-based
						int timing = groupEvents[groupIndex].Time + delay;
						if (timing > maxTiming) maxTiming = timing;

						VirtualEvent vevent = VirtualEvent(
							groupEvents[groupIndex].NeuronIndex,    // presynaptic neuron
							connection->PostNeuron,                  // postsynaptic neuron
							delay,
							connection->Weight,
							1                                       // link level
							);

						virtualEvents[timing, numPostspikes[timing]++] = vevent;
					}
				}
			}
		}

		// attempt to extend
		int t_last = maxDelay + maxDelay + latency + 1;
		for (int t = 0; t < t_last; t++)    // note that t_last can be increased within the loop
		{
			while (numPostspikes[t] > 0)
			{
				numPostspikes[t]--;
				int currentPostNeuron = virtualEvents[t, numPostspikes[t]].PostNeuron;

				// check if already processed
				bool alreadyProcessed = false;
				for (int j = 0; j < numFired; j++)
				{
					if ((currentPostNeuron == groupEvents[j].NeuronIndex) && (groupEvents[j].Time > t - 20))
						alreadyProcessed = true;
				}

				// check inhibition list (excitatory only)
				if (!alreadyProcessed && (currentPostNeuron < numExcitatoryNeurons))
				{
					for (int j = 0; j < numInhibitionList; j++)
					{
						if ((inhibitionList[j].NeuronIndex == currentPostNeuron) && (inhibitionList[j].Time <= t) && (inhibitionList[j].Time + inhibitoryConvergenceWidth >= t))
							alreadyProcessed = true;
					}
				}

				if ((alreadyProcessed == false) && (numFired < maxPNGLength))
				{
					int numLinks = 0;
					int excitationCount = 0;
					int startConvergentOverlap = -1;
					int endConvergentOverlap = 0;
					int firingTime = t + 4;

					// check how many virtual events converge on currentPostNeuron around time t
					// within a temporal overlap of temporalConvergenceWidth
					// the actual convergence on currentPostNeuron is in the range startConvergentOverlap to endConvergentOverlap
					for (int tf = 0; tf <= temporalConvergenceWidth; tf++)
					{
						for (int j = 0; j <= numPostspikes[t + tf] - 0.5 * tf; j++)
						{
							VirtualEvent currentEvent = virtualEvents[t + tf, j];
							if (currentEvent.PostNeuron == currentPostNeuron)
							{
								if (startConvergentOverlap < 0)
								{
									startConvergentOverlap = tf;
								}
								endConvergentOverlap = tf;
								convergedEvents[excitationCount++] = currentEvent;
								firingTime = t + tf;
								if (numLinks < currentEvent.LinkLevel)
								{
									numLinks = currentEvent.LinkLevel;
								}
							}
						}
					}

					if (
						((startConvergentOverlap + 1 >= endConvergentOverlap) && (excitationCount >= minConvergentSpikes - slacke))
						|| (((currentPostNeuron < numExcitatoryNeurons) && (excitationCount >= minConvergentSpikes))
						|| ((currentPostNeuron >= numExcitatoryNeurons) && (startConvergentOverlap + 2 >= endConvergentOverlap) && (excitationCount >= minConvergentSpikes - slacki)))
						)
					{
						// add new firing events
						int lastEventIndex = numFired++;
						groupEvents[lastEventIndex].NeuronIndex = currentPostNeuron;
						groupEvents[lastEventIndex].Time = firingTime + latency;
						if (excitationCount == 2)
						{
							groupEvents[lastEventIndex].Time = groupEvents[lastEventIndex].Time + latency;		// longer latencies for weaker inputs
						}

						// add new links
						for (int j = 0; j < excitationCount; j++)
						{
							VirtualEvent convergedEvent = convergedEvents[j];
							LinkedEventPair^ eventPair = gcnew LinkedEventPair(
								// this will only approximately match the presynaptic neuron firing event (depending on jitter)
								FiringEvent(groupEvents[lastEventIndex].Time - convergedEvent.Delay - latency, convergedEvent.PreNeuron, EventLabel::Foreground),
								FiringEvent(groupEvents[lastEventIndex].Time, groupEvents[lastEventIndex].NeuronIndex, EventLabel::Foreground),
								convergedEvent.Weight,
								convergedEvent.Delay,
								numLinks);

							links->Add(eventPair);
						}

						// try to add more virtual events
						Dictionary<int, List<Connection^>^>^ outgoingConnections = groupedConnectionDictionary[groupEvents[lastEventIndex].NeuronIndex];
						if (groupEvents[lastEventIndex].NeuronIndex < numExcitatoryNeurons)
						{
							for each(int delay in outgoingConnections->Keys)
							{
								for each(Connection^ connection in outgoingConnections[delay])
								{
									// if the connection is strong ...
									if (connection->Weight > weightLimit)
									{
										//timing = groupEvents[i].Time + delay + 1;     // in the original, delay is zero-based
										int timing = groupEvents[lastEventIndex].Time + delay;
										if (timing > maxTiming) maxTiming = timing;

										// BUGFIX: it is possible for 'timing' to exceed the maximum index of the _postspikes arrays
										// the next line checks this limit
										if ((timing < maxPNGLength) && (numPostspikes[timing] < totalNeurons))	// check the limit
										{
											// tentatively add to the group
											VirtualEvent vevent = VirtualEvent(
												groupEvents[lastEventIndex].NeuronIndex,    // presynaptic neuron
												connection->PostNeuron,                      // postsynaptic neuron
												delay,
												connection->Weight,
												numLinks + 1                                      // link level
												);

											virtualEvents[timing, numPostspikes[timing]++] = vevent;
										}
									}
								}
							}
						}
						else
						{
							for each(int delay in outgoingConnections->Keys)
							{
								for each(Connection^ connection in outgoingConnections[delay])
								{
									// in the original, delay is zero-based therefore subtract 1 from time
									inhibitionList[numInhibitionList++] =
										FiringEvent(groupEvents[lastEventIndex].Time + delay - 1, connection->PostNeuron, EventLabel::Foreground);
								}
							}

							// find old items on the inhibition list
							int j = 0;
							while (inhibitionList[j].Time + inhibitoryConvergenceWidth < t) j++;

							if (j > 10)	// needs cleaning
							{
								for (int k = j; k < numInhibitionList; k++)
								{
									inhibitionList[k - j] = inhibitionList[k];
								}
								numInhibitionList = numInhibitionList - j;
							}
						}

						if (numLinks > maxLinkLevel)
							maxLinkLevel = numLinks;

						if (t_last < maxTiming + 1)
						{
							t_last = maxTiming + 1;
							if (t_last > maxPNGLength - maxDelay - latency - temporalConvergenceWidth - 1)
								t_last = maxPNGLength - maxDelay - latency - temporalConvergenceWidth - 1;
						}
					}
				}
			}
		}

		if (t_last == maxPNGLength - maxDelay)
		{
			for (int d = t_last; d < maxPNGLength; d++)
				numPostspikes[d] = 0;
		}

		//Console.WriteLine("{0} {1}", NL_max, t_last);
		if (maxLinkLevel >= minPNGPathLinks)
		{
			// Check that each neuron in the triggering triplet has links to at least two excitatory postsynaptic neurons
			array<int>^ used = gcnew array<int>(groupSize);
			bool discard = false;
			for (int i = 0; i < groupSize; i++)
			{
				used[i] = 0;
				for (int j = 0; j < links->Count; j++)
				{
					if ((links[j]->PreEvent.NeuronIndex == groupEvents[i].NeuronIndex) && (links[j]->PostEvent.NeuronIndex < numExcitatoryNeurons))
						used[i]++;
				}

				if (used[i] == 1)
					discard = true;
			}

			if (discard == false)
			{
				numPolychronousGroups++;

				if (verboseLevel >= 1)
				{
					Console::WriteLine("{0} {1} {2} {3} {4} {5}", tag, numPolychronousGroups, anchorNeuron, numFired, maxLinkLevel, groupEvents[numFired - 1].Time);
				}

				// throw a 'found' event or write to the console
				this->OnFoundPNG(tag, numPolychronousGroups, anchorNeuron, numFired, maxLinkLevel, groupEvents[numFired - 1].Time, groupEvents, links);
			}
		}

		return numPolychronousGroups;
	}
}
